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ABSTRACT 



The periodic position and velocity perturbations of an 
artificial earth satellite are developed to the first order for all 
J , based on the theory by Brouwer as extended by Giacaglia. An 
explicit formulation is also provided for the subset J^ , J„, J, , 
The use of a position and velocity formulation circumvents the 
equatorial and circular orbit singularities found in conventional 
developments. The definition of the mean elements of the theory 
is modified to reduce the complexity of the position perturbations, 
as suggested by Merson's Theory, and the resulting changes to the 
secular terms are developed. In order to facilitate an empirical 
correction for drag, the observed mean motion is introduced as a 
mean element in place of the semi -major axis. 
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SECTION I 
INTRODUCTION 

The motion of a near-earth satellite is, in the first approximation, 
"Keplerian"; i.e., it conforms to certain empirical laws formulated in 
the seventeenth century by Kepler. In his Principia , Newton demonstrated 
that these laws described motion in an inverse-square force field. The 
force (or negative potential) function for such a field is of the form 

u=^ (1) 

r ^ 

where r is the geocentric distance of the satellite and p, is the 

gravitational constant. The motion is conventionally described by six 

"orbital elements," a, e,, I, M, uu , Q. Kepler's second law states that 

the motion occurs along an ellipse with one focus at the primary. The 

inclination, I, and the argument of the ascending node, Q, serve to 

locate the plane containing this ellipse. The eccentricity, e, and the 

argument of perigee, UU , define the shape of the ellipse and its 

orientation within the orbital plane. The semi-major axis, a, provides 

the scale of the ellipse as well as the orbital period; from Kepler's 

third law the period P is given by 

P . ^ a''" (2) 

The location of the satellite within the ellipse is given by the mean 
anomaly, M, which measures the area swept out by the radius vector 
since perigee passage. In accordance with Kepler's first law the area 



swept out and hence the mean anomaly increases at a uniform rate; then 



M(t) = M + n t 
o o 



where n is the mean motion, given by 



(3) 



n = 
o P 



2t\ _ h -3/2 



M^ 



(4) 



The mean anomaly must be converted to a geometric angj.e to be of use; 
the eccentric and true anomal ies , E and v , are related to M by 
Kepler * s equation 

E - e sin E = M (5) 

which must be solved by iteration, and by 



V 1-fe 



E 
tan - 



Following these computations, the geocentric position and velocity 
vectors _r and _r are given by (See Figure 1) , 

r = r U 

_r = r U + rvV 

where 

r = 



r = 



rv = 



1+e cos V 






^' 


e sin 


V 




?.' 


(1+e 


COS 


v) 



(6) 



(7) 
(8) 

(9) 

(10) 

(11) 
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COS u COS n - sin u sin fi cos I 
cos u sin n + sin u cos fi cos I 
sin u sin I 

- sin u cos Q - cos u sin Cl cos I 

- sin u sin Cl + cos u cos cos I 
cos u sin I 



(12) 



(13) 



with 



P = aTl 

r\ = (1-e ) 



u = v-Hu 



(14) 
(15) 
(16) 



This simplified model cannot adequately represent the motion of a 
satellite; it must be augmented to include the effects of various 
"perturbing forces/' These forces may be conservative, i.e. , gravita- 
tional, or may affect the energy of the satellite's orbit, e.g., atmo- 
spheric drag and solar radiation pressure. More sophisticated models 
may be implemented via numeric integration or *'specia!i. perturbations" 
techniques, which fall into three major categories: 

(a) Integration in cartesian coordinates of accelerations, 
resulting from all forces acting on satellite^ to obtain 
position and velocity (Cowell's method). 

(b) Integration in cartesian coordinates of accelerations 
resulting from perturbing forces only, to obtain deviations 
in position and velocity from a Keplerian orbit (Encke's 
method) , 



(c) Integration of the variations due to perturbing forces in 
the orbital elements of an "osculating" Keplerian orbit, 
i.e., the Keplerian orbit defined by the position and 
velocity of a satellite at each instant (Variation of 
Parameters) . 

These methods have a theoretical disadvantage, in that the 
accumulation of roundoff and truncation errors must eventually result 
in inadequate precision; it appears that in practice the length of arc 
is limited by the uncertainty in modeling non-conservative forces, which 
involve not only a complex and highly variable atmospheric structure, 
but also the configuration and orientation of the satellite. A more 
practical problem in some applications is that the integration must be 
carried from epoch to the most distant observation, regardless of 
whether useful data exists in the intervening period. In addition, the 
integrated orbit , whether in terms of coordinates or elements , provides 
little insight into the effects of the various forces operating on the 
satellite, so that it is difficult to identify and correct deficiencies 
in the model. Despite these disadvantages , special perturbations 
programs are widely employed for precision tracking where the frequency 
of data mitigates their relative inefficiency, or the cost is justified 
by the requirements for maximum precision. They are also extensively 
employed in feasibility studies and similar investigations, to avoid 
time consuming (and possibly impractical) analytic developments. 



It. is generally possible to obtain analytic expressions for the 
effects of the perturbing forces, to any desired precision and for any 
time span. Such "general perturbations" models are universally employed 
in routine cataloging systems, where a considerable number of satellites 
must be tracked with data that is sparsely distributed in time. In 
addition, general perturbations are usually employed in satellite 
geodesy. A considerable number of analytic theories have been developed 
for the conservative perturbing forces, i,e. , the departure of the 
earth's gravity field from an inverse-square law and lunar and solar 
tj;ravity- 

For greater efficiency, semi-analytic theories are often employed 
'.or the luni-solar perturbations which are relatively small and of low 
frequency; the same approach is generally followed for the solar 
•adiation pressure perturbations. Thus, for solar radiation pressure, 
..:he analytic development may be carried through a fornal integration of 
■:he perturbations, but the results are left as a function of the limits 
)f integration. These limits depend on the points at which the satellite 
enters and leaves the earth's shadow, which vary slowly with time- The 
valuation of the perturbations proceeds n revolutions at a time, 
A7i th the shadow limits reevaluated at each step. 

Analytic models of the drag perturbation have been produced for 
simplified atmospheric models . The theory is complex, particularly 
when interactions with the earth's oblateness perturbations are considered. 



As a result, empirical models are generally employed, with only the long 
term effects of drag considered. The results are generally satisfactory 
for high altitude objects^ but there appears to be considerable merit 
in the development of a semi-analytic drag theory. 

This paper deals only with perturbations due to the earth*s gravity. 
In geodesy the gravity field is described in terms of a reference 
ellipsoid, a reasonably tractable figure which approximates the figure 
of a rotating fluid in equilibrium to about 1 part in 10 , The actual 
gravity at any point is shown in terms of a map of the elevation or 
depression of the *'geoid*' with respect to this ellipsoid; this "geoid** 
is an equipotential surface, i.e., a surface everywhere perpendicular 
to the local vertical. Before artificial satellites were launched the 
ellipsoid and "geoid'* were determined from the reduction of direct 
gravity measurements and from astronomical determinations of the 
deviation of the local vertical from the local perpendicular to the 
ellipsoid. This "geoid** data is not employed in the theory of an 
artificial satellite, however. An analytic expansion for the potential 
is required; in spherical polar coordinates the generalized force (or 
negative potential) function is a series of Legendre polynomials and 
associated functions: 



Z I— J F (sin 3) + S J F (sin p) cos m iX'\ ) 
- I r ( n n . nm nm nm 

n=2 \ Tii=l 

n 
0° n a 



1 + E E {-zrl P_„ (sin 0) (C^„cos niK + S^^sin m\) 

(17) 



r ) „ , \ r f nm nm nm 

n=2 m-O 



where a is the equatorial radius and 6 , \ are the geocentric 
latitude and longitude of the satellite. The Legendre polynomials and 
associated functions are defined by 



. l.mll jn+m . l.m/l I(— 7-) , ^\^r^ n.\i n-m-2j 

fy\ = (1-^ ) d 2 ^ n ^ (1-x ) ^ (-l)-'(2n-2j): x 

2^: dx""^ 2" j=0 J- ("-J)- (n-m-2j):" 



(18) 

where I(— r— ) is the integer part of {'^TT') * The P (or P ) harmonics 
2. Z n no 

are "zonal/' while the P harmonics are "tesseral," The largest 

nm 

-3 
coefficient is J ; it is of order 10 . The remaining coefficients do 

not exceed the second order, i.e,, 10 , In order to evaluate the effects 

of the harmonics, it is necessary to substitute orbital parameters for 

r, 3, and \. In general, a method of successive approximations must 

be employed, so that a series of perturbations of increasingly higher 

order arises, e.g., 

first order: J 

2 
second order*: Jn 3 J j ^ 

2 n nm 

3 
third order: J- , J-J , J^J 

2 2 n 2 nm 

Most general perturbation theories neglect periodic effects of the 

second order; the residual perturbations will then be on the order of 

15 meters. However, under certain circumstances the perturbations due 

to higher order terms are amplified and must be included in a first 

order theory. If the potential function and its derivatives are 

expressed in terms of conventional orbital parameters, they will be 

found to have arguments of the form: 



COS 



sin 



(n-2p+q)M + (n-2p>JL) + m(Cl-X) 



where n,m are the indices of the harmonic, p ranges from to n 
(it is the parameter of a power series in sin I and cos I) , and q 
ranges from -oo to «» (it is the parameter of a power series in e, 
with the lowest power of e being e' ' ) . From a simple first order 
theory, it will be found that M, uu , and Q all increase linearly 
with time, so that when the perturbations are integrated divisors will 
arise of the form 

(n-2p+q)n + (n-2p)uu + m(n-X) 
o 

where n is the perturbed mean motion. The perturbations are classi- 
o 

fied in terms of n, m, p, q as follows: 

(a) Secular terms p = (n+q) /2 

q = 
m = 

These terms give rise to a linear increase in the elements 
H, uu , fi, and are therefore computed to the second order 
in a first order theory (so that the theory is valid for 
about 10-20 days, after which the neglected third order 
terms exceed the second order) , These terms only arise for 
even order zonal harmonics, i.e., n=2, 4, ... ; the 
values of the even zonal harmonics are generally based upon 
observed secular perturbations . 



(b) Long period terms p = (n+q) /2 

q 9^ 
m = 

These terms have a divisor of the form -quo which is of 

-3 
order 10 ; second order forces therefore integrate into 

first order perturbations and must be included in a first 

order theory. There is no J term of this form; if there 

were, a different type of solution would be required (there 

2 
is a J^ term of this form which reduces to J on 

integration). There is a special case for the "critical 

inclination" I ?« 63.4 , where UJ is of order 10 , so 

that a "resonance" occurs. In this case, either a special 

solution is employed or the long period terms are not 

integrated, i.e., they are left in the form of secular rates. 

The long period perturbations for even zonal harmonics are 

factored by the eccentricity e and can often be ignored; 

this is not the case for the odd zonals whose values are 

usually determined by analysis of observed long period 

variations in eccentricity and inclination. 

(c) Short period terms p ^ (n+q) /2 

These terms have a divisor containing n so that the order 
of the perturbation remains unchanged upon integration. 
Therefore, only the J terms need be included in a first 
order theory. 
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(d) Tesseral harmonic terms m ^ 

There are two cases of interest here. For p = (n4-q)/2 

there are terms with frequencies near some multiple of the 

siderial rate, since to the zeroth order 

(n-2p)LiJ + m(n-X) ^ -mX. 

The integration results in an increase on the order of 

n /mX or about 16/m for near earth satellites. These 
o 

terms contribute perturbations on the order of 100 meters, 
and decrease in importance as n and m increase . 

For 

(n-2p+q)n 4- (n-2p)uu 4- m(n-X.) ^ 

there is a resonance analagous to that holding near the 

critical inclination. The resonance will be in general 

larger for smaller values of (n-2p) and q. The principle 

resonances thus arise for 

m «;^ n n odd 

o 

where n is expressed in revolutions /day . Obvious cases 
o 

of potential near resonance are 24 and 12 hour satellites. 
High order resonances , e.g., (n,m) of (13 , 13) , (15 , 13) , 
and (15,14) have been reported for certain satellites with 
magnitudes on the order of 100-150 meters and periods of 
2.5-5 days . Obviously, by going to a sufficiently high 
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order harmonic, a resonance can be found for any near 

earth satellite. Fortunately, the net effect of these 

higher order terms is considered to approach the second 
order. 

This paper does not deal with the tesseral harmonic perturbations. 
It is limited to some minor modifications of the secular terms developed 
by Brouwer and extended by Giacaglia , and to a non-singular 
development of the long and short period perturbations due to the zonal 
harmonics. In the Brouwer and Giacaglia papers the perturbations of 

conventional elements are computed, which leads to singularities for 

(3) 

low eccentricity or inclination Lyddane showed that the problem 

could be circumvented by either computing the perturbations to "non- 
singular" elements, e.g., e cos M and e sin M, or by computing the 

perturbations in the position and velocity vectors. The former approach 

(4) 
is employed in most general perturbations ephemeris generators , while 

the latter approach is employed in this paper. Although Garfinkel , 

Kozai , and Merson have computed some of the position perturbations, 

velocity perturbations have generally been neglected. 

The use of position and velocity perturbations has the advantage 
of revealing the "real" or observable effects of the perturbing forces; 
Merson , for example, has shown that some of the apparent perturbations 
of the orbital plane affect only the velocity vector and can be ignored 
in a tracking network based on positional data. In addition, a 
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position and velocity theory appears to be somewhat more efficient than 
a "non-singular" elements theory, particularly when only positional 
data is used for element correction. The position and velocity theory 
has one disadvantage, in that the frequency of the long period terms 
becomes comparable to the short period terms, and they must be recom- 
puted for each ephemeris point. (However, they are recomputed for 
each point in most theories, whether or not the computation is necessary.) 
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SECTION II 

PERTURBATIONS IN POSITION AND VELOCITY 



In developing the perturbations it is convenient to use the angular 
momentum unit vector W, given by 



W = U X V = sin n sin I 

- cos n sin I 
cos I 

The perturbed position and velocity may be computed as 

r = (r + 6r) (U + 6U) 

r = (r + 6r) (U + 6U) + (rv + 6rv) (V + 6V) 
or the perturbations alone may be computed as 

6r = 6r U + r 6U 

6r = 6r U + 6rv V + r 6U + rv 6V 
ignoring second order terms. 

The quantities 6U and 6V may be written as 

6U = V (6u + cop I 60) + W (sin u 61 - cos u sin I 60) 

6V =-U(6u + cos I 60) + W (cos u 61 + sin u sin I 60) 
and hence we have 

6r = 6r U + r (6u + cos I 60) V 

+ r (sin u 61 - cos u sin I 60) W 



6£ = 6r - rv (6u 
+16 rv + r (6u 
+ r (sin u 61 
+ rv (cos u 61 



+ c 



cos I 60)1 U 
OS I 60)1 V 
cos u sin I 60) 
sin u sin I 60) 1 



(19) 



(20) 
(21) 

(22) 
(23) 



(24) 
(25) 

(26) 



W 



(27) 
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For Or, 6r^ and 6u we can either use Taylor series expansions 
in the conventional elements , or the ingenious equations of Izsak with 
the Brouwer determining function S written in Hill's canonical 
riables {r, G, H | r, u, Q] , The equations are 



va 



c 6a c . a sin v -^^ 

or = r — - a cos v Oe H e OM (28) 

a r| ^ ^ 

A4. r oa .M.. (1+e cos v) . c . /^^^ (1+e cos v) -u. 

or = - - — + (") ^^ :: ^— sin v 6e + (-) -^» ^— cos v e 5M 

r\ T] 

' (29) 

2 
(2+e cos v) . - , (1+e cos v) -^^ ^ ... ,-^. 

6u = -^' 2 ^^^ V Oe + -^ r ^— OM + OOJ (3o) 

ri r\ 

or 

6r = -^ (31) 

^ ar (32) 

^" = -9^ (33) 

The second set of equations appears much simpler, and has been 
solved for the short period terms by Izsak. They are not so simple , 
however, when dealing with the long-period terms containing tifigonometric 
functions of 'JJ . 

The perturbation in rv may be computed from 



o rv = rv 



1 5a , cos v (1+e cos v) - e ^ sin v (1+e cos v) c^, 
— + i-— ^ 6e - ^^5 e CM 

2 a 2 3 

Tl T] 

(34) 



Brouwer and Giacaglia employ the Delauney canonical variables 



[L, G, H 1, g, h], where 



L = (na) ' 


1 = M 


G = LT] 


g = 10 


H = G cos I 


h =n 
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or 



. . .,G. . ,6G 6r. 
6rv =6(7) = rv (— -) 



using intermediate results for 6G. 



(35) 



If the Taylor series expansions are employed, it is possible to 
rewrite Equation (27) using Equations (29), (30), and (34) as 



e sin V 6a (1+e cos v) sin v ^ 



- (1+e cos v) g^ _ ^^^^ ^^^ ^^ ^^^ _!_ ^^g I 60) U 
ri ' 

( 1+e cos V 6a ^ cos v+e . ^ . ,.,„ ^ ^ .r^. \ .. 

- h ^ 6e + e sm v (6uu + cos I 5w) V 



(cos u + e cos (jj) 61 + (sin u + e sin uu) sin I 60 W 

(36) 



16 



SECTION III 



SHORT PERIOD PERTURBATIONS 



Izsak has already computed 6r, fir, and 6u as 



6r = -7-^ <(sin^ I cos 2u + (1-39^) 



1 + 



21L 



e cos V 



1+e cos V 1-Hl 



(37) 



h ^2% 2 2 

6r = - |Ji r-rr / 2 sin I (l+e cos v) sin 2u 

4 p I 

-h (1-302) ^ ^.^^ j^^ (14-e^os v)2 j I 

J a M 2 

2 ^ "^ 6(1-50^) (v-M) + 4(1-60^ + \'2^ ) e sin v 



(38) 



6u = 



8 p- 



1+ri 



+ (1-30^) (1-ri) sin 2v + 2 (50^ - 2) e sin (2u-v) 



+ (70^ - 1) sin 2u + 20^ e sin (2u-hv) > (39) 



where 



J = coefficient of the second zonal harmonic 

a = earth's equatorial radius 
e 

= cos I 



From Brouwer*s theory, with 

J.a^ 
Y =-2-^ 
2 2a2 
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we have 



— = 3 TT sin I / cos 2u + e cos (2u-v) + — cos (2u-h/) ) (40) 

4 p 

^2% 2 
cos I 60 = - — 0^ \6 (v-M + e sin v) - 3 sin 2u 

4 p^ 

- 3 e sin (2u-v) - e sin (2u+v) \ (4l) 

f 

61 = TT- sin 19(3 cos 2u + 3 e cos (2u-v) + e cos (2u-K7) 

4 p^ 

(42) 

Hence 

c . -x 9.rs 2 e // . 2 . sin 2u ^ • /o \ 

6u + cos I 6n = — <(sin I) 1 — + 2 e sm (2u-v) 

4 p I \ 

- (1-30^) 3 (v-M) + 2 e sin v (|^) 

+ ^ sin 2v I > (43) 

sin u 6l - cos u sin I 6Q = — sin 10 V 

4 / 

- 3 sin u - 4 e sin 'ju 

+ 4 e cos u sin v + 6 cos u (v-M) ) (44) 

J.a ^ I 

:os u 61 + sin u sin I 60 = 1- sin 10/3 cos u + 4g cos (u+v) 

A P^ [ 

- 6 sin u (v-M) \ (45) 
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6 rv = M- 1/2 (1+e cos v) V 



4 P 
(sin^I) 



2 cos 2u + 2e cos (2u-v) + e cos 2u cos v 



(1-30^) j (1+Ti) + e cos V (j~p +^ cos 2v 



(46) 



Substituting these terms in Equations (26) and (27) gives, after 
some simplification, the short period terms in Table I. 

In some tracking programs the relationship of the mean semi -major 
axis and the secularly perturbed mean motion is taken from Kozai*s 
equation 14 



T^^ - 



Since 



n = n 



^ 1 + 3 -^-y- (1-39^) r] 
\ 4 p 



■^2^e 2 
1 - 3 ^ (1-39 ) r\ + .. 

2 P 



(47) 



(48) 



this implies 



a = a 



1 + 3 I (1-39^) r\ + ... 
4 p 



(49) 



where a is the mean semi-major axis of the Brouwer theory, defined by 

n^ a-^ = M- (50) 

o 
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The use of a in place of a in the computation of unperturbed _r 
and £ requires that a compensating perturbation be applied to &£ 
and (}V: 

-~ = - 3^ (l-30') ^ (51) 

a a ,2. 
4 p 

This results in the following changes in Table I: 

(a) The term 2x\ in the U component of 6^ 
becomes - T] 

(b) The term -r] e sin v in the U component of 6r 

r\ 

becomes 4- — e sin v 

Merson has developed a theory in which the short period position 
perturbations are minimized. His formulae (153-157) relate the osculating 
elements to conditions at the ascending node , to the second order in J^ . 
By eliminating all first order terms whose argument is a multiple of uu , 
and terms factored by u , a set of first order pseudo short period terms 
is obtained. These result in the position perturbations given in 
Table II. (The first order terms factored by u in Merson^s theory are 
actually the sum of secular terms factored by M and short period terms 
factored by (u-M) . They can therefore be included in the secular term 
computation, and we have taken this approach.) 

Now, if we define ^, to be Brouwer^s mean elements updated for 

I 
secular perturbations, and €. to be "smoothed" mean elements updated 
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for secular terms, then; 



rCe,) + 6r„ = r(e;^-) + ) -^ Cr *-r ) + 6r 
-^ ^^ -Brouwer -^'' i^ /_^ g^ ve^ r.^; -»- ^r. 



Smoothed (52) 



i=l 



so that we have three simultaneous equations: 



/ .he. 1 -1 



i=l 



6 r ~ ^I. 

Brouwer Smoothed 



(53) 



where 



^e , = e ! - e . 
1 1 1 



Using Equations (26), (28), (30) and Tables I and II, 



Aa A _ . a s in v ^;,^^| _ „ ^ /i -^n^xV 



r — - a cos V Ae + 
a 



eAM = r a, (1-39 ) 
T) / 1 



l+3n+e 2-Hn 

2 1-Kl 



e cos V + — cos 2v 



((2 sin v+l sin 2v) ^ + 11±^^°5-^ Am 4- /^ + 9An)= - ra^(l-392))( 



3 (v-M) + 



2-m 



l-4r| 



2e sin v + — sin 2v 



r (sin u Al - cos u sin lAfi ) = r a sin 19 



6 cos u (v-M) - 3 sin u 



(54) 
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If the Ae are restricted to contain only constants or terms with 
(v-M) as angular arguments these equations may be solved to yield: 



^ = 2 a, (1-30^) r\ 
a i 



Ac = - a^ (1-30 ) 



2+r] 



1-Hl 



er\ 



AI = - 3 a sin 10 



AM = 



Ao) = - 3 a^ (1-50 ) (v-M) 



An = - 6 a^ (v-M) 



(55) 



The resulting changes in the velocity perturbations may be obtained 
from an equation similar to Equation (53); the signs of the Ae . 
should be changed to yield: 



6r - 6r 

Smoothed Brouwe 



i=l 



(56) 



With the aid of Equation (36) 



6r - 61 

Smoothed Brouwer 



'^1' a, X 



(1-39") 



e sin V -r| + 



2+Q 
1-K1 



(1+e cos v) 



+ (1-30 ) 



2 f„ . 2+11 



e cos V T| + 



+ 3(l+e cos v)(v-M) 
3e sin v j(v-M) 



i+n 



- 3 sin 10 



cos u+e cos UU + 2(sin u+e sin uu) (v-M) 



U 

y 

W 
(57) 
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which, when added to 6_r from Table 1, yields Lhc velocity perturbations 
in Table II. 

These position and velocity perturbations represent a useful 

simplification of the equations in Table I; corresponding modifications 

2 

to the J^ and J secular terms are derived in Section V, 
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SECTION IV 
LONG PERIOD PERTURBATIONS 



The first order long period perturbations for J , n> 2, have 

(2) (9) 

been developed by Giacaglia and Garfinkel and McAllister 

Giacaglia's results have been employed since they are more readily 

truncated; the published paper unfortunately contains a number of errors: 

(a) For the even harmonics, the upper index limit for k in 

2 is incorrect, since values of i within the permitted 

kij 

P 
range do not exist for k= ■^. Furthermore, for k = 0, i 

can exceed j; in this case, 6 ,, ^, ^. vanishes.'^ 

' p+1, 2j, 2i 

Therefore 

2 mm Y^ — > 2 
kij k=0 i=I j=l 

Similarly for the odd harmonics 



2 Mn^-E-p^) 2 
kij k=0 i=0 j=0 



p-2i -l . , . 
where -^ — r-' — was incorrectly given as the upper index 



limit of i. 



These practical proljleius reflect theoretical errors; i.e., ^^ ^ 2 *^^^^'^^ 
nrc secular terias aid also appear in the >: portion of the disturbing 
function. kj 
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(b) In equation 21, the term 



2 
e 

should be 

e 

(c) Equation 22 must be multiplied by H, 

(d) In equation 28, the factor 

2i+l 
should be 

2i+l 

(e) In equation 30 the term 

sm I 

should have a minus sign, and the factor 

cos (2i+l)R 
21+1 

has been omitted. 

Also, Giacaglia has the wrong signs for the third through fifth 

harmonic Brouwer coefficients (k , v > ^"^^ A ^) as functions of J . 

^ n* ^ n n.O n 

With these corrections, Giacaglia's results may be summarized as 

6e = S Ae (58) 

where S and ^e for even and odd n are summarized in Table III. 
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It is convenient to reduce these equations to a common form. First , 
note that the summations can be written as 

n-2 n-2 n-2i 

2 2 2 

H '\; Z n even 

i=l j='i k=0 

n-3 n-3 n-2i-l 

2 2 2 

Z 2 S n odd 

i=0 j=i k=0 

Now introduce a variable \ 

X = 2i n even 

= 2i+l n odd (59) 

so that 

X = 2, 4, ... , (n-2) n even 
=1,3, ... , (n-2) n odd 
Similarly introduce 

li = 2j n even 

= 2j+l n odd (60) 

so that 

^, = X, X+2, ... , (n-2) 
Also introduce 

V = n-2k 
with 

V = X, ?.+2, . . . , n (61) 



Note that S is factored by e sin I. Moving this factor from S 
into the Ac , and employing 

% = n/2 - lu (62) 

we obtain the results in Table IV. 

It is now a simple matter to obtain 



6r = -aTj sin I S -j cos v co 



s \^ + -- sin V sin \ 



5} 



(63) 



f^ • M* 1/2 2 f p* 1 

^r = (— ) (1+e cos v) sin I S -j sin v cos ^?-r cos v sin X? > (64) 



6rv = (— ) (1+e cos v) sin I S 



\/ \ (1+e cos v) cos V cos X^ "*" T" sin v sin X§ - e cos X§ > 



(65) 



6u + cos I 6n = sin I S ^ (2 sin v+e — ~ ) cos X^ 



ro j_ cos 2v . 3e . M' sin \g 
- (2 cos v+e — 2 ~ ^ X 



+ e (2n-5) ^^^ } 



(66) 



sin u 6 I - cos u sin I 6 fi = -e S 



X, \ sin u cos XZ + cos u | v - 



:os u 6 I + sin u sin I 6 Q = -e 9 S 



2 \ 
10 sin I 

2 
1-59 ' 



sin 



XZ 



(67) 



X \ cos u COS X5 "*" sin u 



^ 10 sin^I 



- V 



sin \g 



\ 1-50 

from which, with Equations (26) and (27), the results in 
Table V are obtained. 



f^} 



(68) 
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It will be noted that S vanishes for n=2. If this \;ore not the 
case, it v;oulcl be necessary to modify the basic unperturbed solution, 
since the resulting terms would have a J^/J coefficient and could not 
be treated as a perturbation. However, the maj^nitude of J is such 
that the long period solution for J^ must be carried to one higher 

approximation than for the other J , yielding perturbations with the 

2 
coefficient J^/ J^ or J^ . As given by Brouwer, these terms are; 



6a = 



(69) 



2 2 
6e = S T] sin I cos 2a) (70) 

6l = -S^ e 9 cos 2'ju (71) 

2 ^ 

6M = 82^ sin I sin 2(jj (72) 



e 



2 
6tu = - S^ e sin I 



1 . 1 9 109 1 

^ e sin I (1-59 ) (1-159 ) J 



2uu (73) 



6Q = - S^ e sin 19 — ^- + ^ r- sin 2ju (74) 

^ I sin I (1-59 ) (1-159 ) / 



where 

J a ^ 
Sl = i- -^ ^-Sln^ (1-159^) (75) 

2 16 p2 (1-592) 

The additional perturbations are easily obtained: 
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6r, = r -77 — — ^ — -y~ (1-159^) ( - sin I (1+g cog v) cos (2u-v) U 

TO • /-o \ J c sin 2u 1 
+ sm I 2 sin (Zu-v) -I r V 



„ . , „ s 10 sin I cos u sin 2ilJ , „ 
e9 sin (u-2UU) ;:; ^i;^ | W 



(1-59^) (1-159^) 



(76) 



^12 = 



M^ 1 '2 % e sin I (^.^3^2^ 



V\ '' p2 (1-592) 



X 



sin I (1+e cosv) sin (2u-v) - 



e sin 2ai 



+ sin I cos (2u-v) + e cos 2iJU 



e sin V sin 2'iJ 



^ ^ rL \ . ., . 10 sin I (sin u+e sin uu) sin 2uu 
e9 cos (u-2uu) + e cos '^ -I- ^^ — ' :r^ 



(1-59^) (1-150^) 



W 



(77) 



In most general perturbations formulations it is customary to employ 

2 
only the J^, J^, and J long period perturbations. For comparative 

purposes, these are given explicitly in Table VI. An evaluation yields 



For J-: 



\ = |j, = 1 V = 1,3 



1 J^ a „ -1 

s-4 vF(i-=^'> X 



4 2 ^ = 1 
•5 sin I V = 3 



(78) 
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«^3 



r ( sin I (1+e cos v) sin u U 

2 ^ 



+ sin I (2-l-G cos v) cos u V 
+ eO cos V W 



(79) 



^^^3=1?) Z37F 



sir I (1+e cos v) cos u U 



- sin I (sin u+e sin (u) V 
- e9 sin v W 



(80) 



For J, 



X = ^JL = 2 V = 2,4 



c = . 2_ 4 e e sm I 

' " j^ ■>' »-5^'' 



b V = 2 

X \-7 sin^ I V = 4 



(81) 



5r = r --7 ^ ^ ^ ^^" •*• (1-79^) < - sin I (1+e cos v) cos (2u-v) U 
J^ P (1-59^) 



+ sin I 2 sin (2u-v) + 



e sin 2u 



- e9 sin (u-2iii) 



2 sin I cos u sin 2ju 
(1-50^) (1-79^) 



W 



(82) 
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'i/, " 



1p) " J, P^ (1-59^) " ' 



j/ i -sin I (1+e cos v) sin (2u-v) 



e sin 2a) 



+ sin I cos (2u-v) + e cos 2tAJ - 



e sin V sin 2uu 



r, / o/iA . . 2 sin I (sin u+e sin UU) sin 
e9 cos (u-2tAJ) + e cos UU + ^ — — ^ 

\ (1-59^) (1-79"^) 



From which the equations in Table VI are easily obtained. 



2UU 



W 



(83) 
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SECTION V 
SECULAR TERMS 



The results of Brouwer and Giacaglia^ "^ are given in Table VII, 



where 



n-2 
. J^a^^^' ^ (-l)N2n-2k): e^^ sin"-^^ 



' 2_ 02" o" ^ ^ ,2(j. 



n*2 
n even 



k-0 j=0 



^^k:(n-k):(^:)^j:)^n.2j-i): 



(84) 



In most practical applications it is the parameter n rather than 
a which is determined. From the relationship 



Y5^ = 1 - a +a + ..• 

where (1-Kl) represents the terms factored by n in the equation for 
n we obtain to the second order: 



n = n (1-a + a + . . .) 

o 

It (1*1 i!|!^ (1.3,2,* ^Wl^e^ 



3-309^ + 359*^ 



t2 ^ 
3 i2!e _ 

"T28~7r^ 



10-&r\-25e^ + i-60+U8r\+90e^)9^ + (130-72n-25e2)9^ 



TT (l - STl (2j ^ - 3) + 



j„ term above 



(85) 
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if n J.S substituted icr n in the equations for w and in 



Table "^l'-; 'c ^^ jc.\,,-i- ...ij of the j" ter-is is reduced; they bee 



ome: 



F^r 



J a ' 
n ilg ---f- |-lO-25e^ + (-36+126e^)0^ + (430-45e^)9^ 



(86) 



For p 



.2 4 



if-9e^ + (-40+5e^)9^J 



(87) 



bj ccuipu'ted from the relationship 



a^-^'^'.-'^' -.l/^n-2/3 a.2/3a- 1/9 a^..) 
o 



1 >L. « . -.n ve i, ju:.^i. 



1/3 ^ 
= U n 



2/3 



( 1 



1/2 -^^ Tl (1-30^) 



0-4 „'t 



10+12n-25e^ + (-6O-72Tl+9Oe^)0^ 
+ (13O+lO8ri-25e^)0'^ 



15 •^4^e ^ 2 

6Z— T^ ^ 



3-300^ + 350^ 



V.^/3n'''Nl+STl(|ii^- 21 + 



J term above 



(88) 
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Note that p, which is a function of a, enters into these equations. 
If we employ 






2 
in the equations in Table VII five J T] terms are added to Equations 

(86), (87), and (88). It is therefore more efficient to compute p 

-2 
from n and p from 



-2 2 

P = P 



^ 2 

^2% 2 I 
1 + "^|"T1 (1-39^) + ... (89) 



^ 2 

for use in the J^ terms. Either p or p may be used in the J 

terms, since the error in using p is of the third order. If, however, 

Equation (88) is only carried to the first order, it may be entered with 

p, and p computed from a, so that Equation (89) is not required. 

If the Kozai a is to be employed as an element, then Equation (88) 
is replaced by 



/i J. 3 2 e 

^ " ^ V ^ 4 T 



r] (1-39^)] 



1/3 --2/3 

\x n 



' J a 2 

1 +-^-|- Tl (1-30^)1 (88a) 

\ 4 p^ 
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t\' the first order. The use of n and 



p = a (1-e ) 



in the equations in Table VII would again lead to additional terms 

in T\ in Equations (86) and (87). The preferred approach is therefore 

to correct "p using 



P = P 



3 ^2^e 
1 +j ' 



1 \ 



(89a) 



A rather more conplicated approach is employed in the GEPERS 

(4") 
ij.r:i..;i developed for use in the SPACETRACK system. The parameters 

ciit^lcyeJ include a and "p, with a mean motion parameter n defined 

t ; 

2 



n^a"* =u ;i +^ ^:r ^ (1-^Q ) 

p 



(90) 



\./. l-:h is interpreted to mean 



n = n 



1 3 -^2% ,. „2, 

1 - -^ — — T\ (1-39 ) 

P 



= n 



1 - \ -^ n (1-39^ 



3 2 e 



^ — T1 (48Tl-288n9 + 432119 ) 
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(where either p or p may be used in terms of order J^-) ^f 

this equation is subtracted from the equation for "n in Table VII 

and n is subtracted for n in the second order terms we obtain: 

o 



t2 ^ 



4 



10-3Zn-25e^ + (-60+192TH-90e^)9^, 
+ (130-288ri-25e^)9'^ 



45 -^4^6 2 

128 -1:4- "^ ^ 
P 



3-309^ + 359'^ 



(91) 



The expressions for ulj and Q in GEPERS employ n and "p in 
place of n" and p; 'n agrees with "n to the first order, but 
the use of 'p makes it necessary to modify Equations (86) and (87) 
The results may be obtained from Equation (86a) and are: 

For uu 



t2 ^' 
_3_ Ve_ 
128 p4 



•10-48ri-25e^ + (-36+384rn-126e^)9^ + (430-72On-45e^)9'^ 



(86a) 



For n 



t2 4 
B9 ^ '^'' 



2s „2 



rr^ / ]4-24n-9e + (-40+72ri+5e )9 

32 M- I 



(87a) 
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Suppose now that we wish to use the smoothed elements n, a', 
e', I', 0)' and Q' as defined in Equation (55). If we enter the J^ 
terms for uj and Q in Table VII with 

p' = a- (l-e-^) 



Tl 



= Cl-e'^)^'^ 



9' = cos I' 



then the J. terms in Equation (86) and (87) will be changed by 

2 



6uu = n 



' 3 ^2% 2 



- 6a ^ 4e 6e , 10 sin 19 61 



r\ 



(1-59^) 



(92) 



6n - n 9 



3 h^e 



2 2 

P 



9 ^3. 4. ^6 Se s in I 61 

a , 2 ■ 

1-e 



(93) 



where the 6e . are related to the Ae . of Equation (55) by 



6e . = e . - e ! = - Ae . 
1 1 i 1 



The results are 



For (JO 



- 3. h\ 
" 128 4 



■42-57e^ + (-20+382e^)9^ + (190-525e^)9^ 



(86b) 
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For Q 



"® 32 T 



■25e^ + (-4+53e^)9^ 



(87b) 



These terms differ from those given by Merson because he employs 

u as the argument of the secular terms rather than M (in the form 

nt )j i.e., his rates apply for a nodal rather than an anomalistic 

2 
frequency and may be obtained by subtracting — and i:- from 

Equations (86b) and (87b). 



Now since from Equation (55), 

2 



J a 

UU' - UU = - 7 -H- (1-59^) (v-M) «-^ (v-M) 
q. z n 



J a * 

n* - n = - 1 -^ e (v-m) «-^ (v-m) 



if we define the increase in mean anomaly to be 



AM = n At 



(94) 



and we define Av to be the related quantity 



Av = AM + (v-M) 



(95) 



(note that Av is not zero at the epoch, unless epoch is defined to 
be at perigee or apogee ) 
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we have 



uu' = uu + uu At +-^ (v-M) 
o n 



(96) 



= uu' +^ Av 
n 



and similarly 



n' = n' +^ Av 

o n 



(97) 



where -r:— and -=:— are obtained from Table VII as modified by 
n n -^ 

Equations (86b) and (87b), It remains to determine a* 
(and, thereby p* ) , Since 

2 



a' = a 



1 +^-^^ (1-39^) 

^ D / 



Equation (88) becomes 
1/3 --2/3 



= M' 



t2 ^ 
^ ^64— T-^ 



10-4X1-256^ + (-60+2411+9062)02 



+ (130-3611-256^)0'^ 



15 Ve_„ 2 
_ —r^ 6 



3-300^ + 350*^ 



1/3 --2/3 



41 lie; 



2 + 



J t6rra abov6 



(88b) 



(n> 2) 
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Thus, to the first order, a' may be computed from Ti using the 
Keplerian relationship and no special provisions are required for 
computing p* . 

The non- conservative perturbations are sometimes represented 
by polynomials in time. If, for example 

M =M + nAt+-^ At + 7 At (98) 

2 6 



li = n + nAt+|At (99) 



then the coefficients for a' 

a' = a + a At + -| At (100) 

are easily obtained from Equation (88b) as 

1/3 -2/3 
a = |jL ' n (101) 

A = - I - a (102) 

3 n 



2 
2 - i 9 In 3 n ■ ^ ^^°^^ 



A similar polynomial can be employed for e* 



i' = e + e At + I At (104) 
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The time dependent terms may be obtained empirically, or from the 
approximation that perigee height remains constant 



li^' "-•> 



= 



which yields 



e = 



1-e 



a = - — 



2 1-e 



(105) 



e 
2 



1-e 
a 



2\ 



a 

;2 



a 
a 



n 



1 n^ 
9 n 



1 .. 

3 ^ 



(1063 



The inclination is usually assumed constant, although it may be 

represented by a similar series. Where these secular perturbations 

on n, a*, and e* exist, they must be reflected in the motion 

of node and perigee. The variations in n are already incorporated 

in Equations (96) and (97). It is usually adequate to carry the 

additional variations due to a', e*, and I* to the first order 

in J^; the necessary derivatives of -=- and — may be obtained 
2 n n 

from Equations (92) and (93). . The resulting equations are 



/ \ 




ua 
o 




n 


n' 


s 


^o 


+ 




^ i 




I > 




[ ) 



Av 1 



r. ^^ 4e 
1-e 



Se + 



10 


sinl9 


1- 


■59^ 


sin I 



SI 



(107) 



where the 6's represent the change from epoch, or 
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= 


f \ 


+ 


n 

4 

n 
^ J 



Av (1+a) 



(108) 



where (using Equations (1023 and (1033 ) 






10 sin 19 
1-59^ 

sin I 



9 



V 



^ (!)^ 



i At 



^ 10 sin 19 
1-59^ 



^ 2 n ^ 2e .. . 

+ -r - + r e + 

3 n ,2 
1-e 



sin I 
9 



\ 



/ 



At 



(109) 



or, for constant perigee height and inclination (see Equations 
(105) and (1063 ), 



3(H-e) n 



'-10-6e in-^ 



l9(l+e) 



!1 ^IRS}-n "' 



In some cases it is desirable to express the secular terms as 
functions of time rather than as functions of Av. It is possible 
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to expand Equations (108), (109) and (109a) using Equation (953, 
and by separating like powers of At to obtain 



I 




M 




'^1 








uu 
2 




/ \ 

uu' 

6 


v"'. 


» 


vPo/ 


+ 


<"y 






+ 


2 
^ J 


At2 + 


a 

6 



At + . 



(110) 



where 



uu 








n 








v". 




a. 

n 






uu 
2 

n 


= 


/ \ 


10 sin 19 




n 




11 . . 4e 
1-e 


1-592 
sin I 


ni 










11 - 5e .' 

6(1+6) "^ 









(111) 



(112) 



(112a) 
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